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Abstract: Laser Speckle Imaging (LSI) is fast, noninvasive technique to 
image particle dynamics in scattering media such as biological tissue. While 
LSI measurements are independent of the overall intensity of the laser 
source, we find that spatial variations in the laser source profile can impact 
measured flow rates. This occurs due to differences in average photon path 
length across the profile, and is of significant concern because all lasers 
have some degree of natural Gaussian profile in addition to artifacts 
potentially caused by projecting optics. Two in vivo measurement are 
performed to show that flow rates differ based on location with respect to 
the beam profile. A quantitative analysis is then done through a speckle 
contrast forward model generated within a coherent Spatial Frequency 
Domain Imaging (cSFDI) formalism. The model predicts remitted speckle 
contrast as a function of spatial frequency, optical properties, and scattering 
dynamics. Comparison with experimental speckle contrast images were 
done using liquid phantoms with known optical properties for three 
common beam shapes. cSFDI is found to accurately predict speckle contrast 
for all beam shapes to within 5% root mean square error. Suggestions for 
improving beam homogeneity are given, including a widening of the natural 
beam Gaussian, proper diffusing glass spreading, and flat top shaping using 
microlens arrays. 
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1. Introduction 

Dynamic Light Scattering (DLS) is a method which utilizes the decorrelation of coherent light 
to estimate the motion of particles in scattering media. DLS is used across scientific 
disciplines to characterize samples in physical, chemical, and biological settings [1-3]. In 
biomedical applications, several DLS methods are used, including: Laser Doppler flowmetry 
[4-6], Diffusing Wave Spectroscopy (DWS) [7-10], and Laser Speckle Imaging (LSI) [11- 
15]. Common to all techniques is that decay of the autocorrelation function of coherent light 
is related to the degree of scattering particle motion. 
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It has been shown that in highly scattering media, the electric-field temporal 
autocorrelation function G\(t) = (E(0)E*(z)) obeys a transport equation [16]. Using this 
model, the autocorrelation function can be calculated and compared to experimental 
measurements in order to determine the dynamics of the light scattering structures. This 
approach has been used with simple, point like sources extensively for fiber based techniques 
such as DWS, in which the distance between the light source and detector is clearly defined. 
In contrast, for imaging based measurements such as LSI, source-detector separations are ill 
defined. However, by taking advantage of the spatial filtering characteristics that turbid media 
enact onto structured light, one can model correlation transport as a function of the spatial 
frequencies in the structured light projected on to the sample. 

This concept has been applied to the transport of light radiance as the foundation for a 
reflectance imaging technique called Spatial Frequency Domain Imaging (SFDI), which has 
been well established and used to determine the absorption, scattering, and fluorescent 
properties in turbid media [17-23]. Recently, we applied experimental methods and modeling 
precedent from SFDI to the transport of correlation by projecting structured coherent light on 
dynamic turbid samples, and measuring the speckle contrast as a function of projected spatial 
frequency [24]. We showed that the speckle contrast follows expected behavior as a function 
of spatial frequency [25], and that correctly modeling the effect using correlation transport 
allows one to fit quantitatively for the dynamical properties of light scattering particles [24]. 

The dependence of speckle contrast on spatial frequency shows that despite the fact LSI 
may be independent of overall intensity of the laser source, spatial variations across the beam 
spot itself can cause the speckle contrast to change. In typical LSI experiments the beam has a 
naturally Gaussian intensity profile, is projected over a finite extent, and may have additional 
artifacts from traveling through a beam spreader such as a lens or glass diffuser. Thus, to 
increase the quality of speckle imaging experiments, it is important to study the nature and 
magnitude of the phenomenon more in-depth. This manuscript first demonstrates effects of 
the beam profile qualitatively using in vivo experiments, then performs quantitative 
predictions and analysis using coherent Spatial Frequency Domain Imaging (cSFDI) of 
controlled liquid phantoms. 

2. Theory 

The speckle contrast parameter is defined as 



where K is the speckle contrast, a is the intensity standard deviation, and (7) is the mean 
intensity. When imaging the speckle pattern with some finite integration time T, K is related 
to the electric field autocorrelation function G\(r) through the expression [14,26] 

* 2 =^f(l-^) 
T J T 

1 0 1 

Here /? represents a multiplicative reduction in contrast not associated with the dynamics of 
the scatterers, such as depolarization, coherence length, or mismatch between the image pixel 
and speckle size. If the scattering objects are static, the field does not change and K reduces to 
p. Assuming the scatterers are in motion, a correlation diffusion equation can be used to find 
G\ analytically as a function of optical properties and the mean square displacement (Ar 2 (r)) 
[10]. However, we and others have found that when imaging speckle over a wide field, the 
diffusion approximation breaks down more readily because of shorter average photon path 
lengths and relatively long camera exposure times [24,27]. Therefore, we utilize Monte Carlo 
(MC) simulation techniques to determine the autocorrelation function. As proposed in 



G,(r) 



0,(0) 



dr. 



(2) 
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[9,10,28], the autocorrelation can be found as a function of separation from a point light 
source, p, by tracking the momentum transfer distribution at discrete p bins and numerically 
integrating over the single scattering decorrelation exponential [29]: 



G 1 (r,p < ) = JP(r)exp 



■fc,?r(Ar 2 (r)) 



dr. 



(3) 



o 



3 



Here P(Y) is the distribution of dimensionless momentum transfer, Y = 1 - cos(6>), and k 0 is 
the wave number inside the medium. 

To take advantage of the cSFDI technique, the autocorrelation is then calculated as a 
function of spatial frequency by applying a numerical Hankel Transform to Gi [24] 



Finally, the speckle contrast modulation transfer function K MTF {j), can be found by inserting 
Gi into Eq. (2). K MTF (f) is a spatial filter that the sample enacts onto the speckle contrast 
image when the coherent source has spatial structure. 

3. Methods 

Unless otherwise noted all experiments were performed using a 30 mw 633 nm HeNe Laser 
(JDS Uniphase, Inc) coherent light source, and images were captured with a Retiga EXI 
cooled 12-bit CCD camera (Qimaging, Inc) at an integration time of T = 15 ms. The optical 
projection axis was angled slightly (10 degrees) to avoid specular reflection. A cross polarizer 
was also placed in the camera lens to help further reject specular light. Speckle contrast was 
calculated using a sliding window mean and standard deviation filter, calculated over small 
sampling neighborhoods of 7x7 pixels. Border areas within 7 pixels of the edge were 
discarded. Speckle flow index was defined as the reciprocal contrast FI = l/K. All 
measurements were carefully designed to meet or exceed Nyquist criteria to minimize 
sampling errors [30]. The value of was found empirically using a static diffusive target and 
calculated to be fj = 0.8. Note that finding /? in this manner may also incorporate shot, readout, 
and other types of noise which can inflate the value. However, with adequate signal and 
cooling this effect should be small, and any error in calculating fj will appear as a 
multiplicative factor in images of speckle contrast and not alter the conclusions of the paper at 
large. 

Two in vivo speckle contrast measurements were performed. In one experiment, images 
were acquired through the skin on a human radial artery, near the superficial palmar branch, 
with a Gaussian illumination profile of FWHM 10 mm. To track location, opaque fiducial 
markers approximately 2 mm in diameter were placed near the artery using adhesive. The 
branch point of the artery was first centered under the peak of the Gaussian, then moved 
diagonally outward towards the tail edge of the beam. A second in vivo measurement was 
taken on a mouse cortex, where Gaussian and flat-top laser beam profiles were used to image 
a 17-month old C57/BI6 male mouse under 1% isoflurane anesthesia. A bilaterial craniectomy 
was performed to expose the brain. 

Experiments were also performed with liquid scattering phantoms of known optical and 
dynamical properties composed of a microbead solution. Phantoms were mixed to a 
concentration of 0.48% w/v with 800 nm Polystyrene beads (Spherotech Inc.) and distilled 
water. The optical properties, calculated using standard Mie Theory and known water spectra 
Motion diffusion characteristics, where the mean square displacement (Ar 2 (r)) = 6D,,r, and D b 
can be calculated using the Stokes-Einstein formula; in this case D h = 6. 1x10" 7 mm 2 /s. 
Speckle contrast data was taken using the following beam shapes: standard Gaussian profile, a 
flat-top profile shaped using a fly's eye lens array pulled from a digital projector unit (Epson, 



G, (*,/) = 2*jG, (t, P )J 0 {lnfp)pdp. 



(4) 



o 
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Inc), and the output through a frosted glass diffusing plate (Thorlabs, Inc). Note that for 
coherent light the intensity pattern formed through a diffuser is speckle in nature and depends 
on the radius of the beam waist and the distance between the diffuser and object [26]. The 
beam width used here was typical for a HeNe output, with 0.5 mm FWHM, and the diffuser 
placed at a distance of 10 cm from the sample. Control measurements were taken with each 
beam shape on a Spectralon reflectance standard. This was done to confirm that variation in 
speckle contrast is not caused by systematic errors such as low light intensity, noise, or 
vignetting. Spectralon is composed of static highly scattering material where no change in 
speckle contrast is expected with respect to beam shape. 

To numerically simulate cSFDI data, a simple time resolved Monte Carlo photon 
propagation algorithm was obtained from the Oregon Medical Laser Center Website [31], and 
a new feature was added to track the weighted dimensionless momentum transfer. Ten million 
photons were propagated through a semi-infinite medium with optical parameters matching 
our liquid phantom above, and P(Y) found at 200 discrete radial bins linearly spaced between 
0 mm and 40 mm. G\(j) was then calculated using Eq. (3) and Eq. (4), and finally the transfer 
function of the contrast image K MTF (f) calculated with Eq. (2). A plot of K MTF (f) is shown in 
Fig. 1 for a typical reduced scattering coefficient in tissue fi/ = 1.7 mm -1 . Note that the 
scattering medium acts as a high pass filter on the speckle contrast. 
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Fig. 1. Plot of speckle contrast as a function of spatial frequency using Monte Carlo 
simulations. This is the "speckle contrast modulation transfer function" K M Tf{f), which is 
applied as a spectral filter on beam profile images to predict remitted speckle contrast. 

The Fourier spectral coefficients of the reflectance image for each beam shape were 
determined by implementing a Discrete Fourier Transform algorithm. For finite exposures, 
the recorded speckle is a stochastic process and thus the mean-filtered image (7) has nonzero 
variance which depends inversely on the number of sampled pixels, determined by window 
size, according to the central limit theorem. This variance manifests itself in the high spatial 
frequencies when using DFT. To decouple the intensity profile's spatial frequencies from this 
natural statistical variance, over 100 independent images of (/) were averaged together, 
effectively increasing the sample number by two orders of magnitude. The diffuse intensity 
image was then represented with a Fourier series, 

(W)^l«„e i(2 * +2 " V) . (5) 

n=0 m=0 
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Fig. 2. a) Human radial artery illuminated with a Gaussian beam, b) Arterial branch placed 
directly beneath the peak (left), then moved diagonally and placed underneath the Gaussian tail 
(right). Dark ellipses represent fiducial markers, c) Visual differences in flow index are seen 
between the two locations. These depend on region of interest, but can reach 13% in the high 
flow arterial sections, d) A profile plot taken horizontally across the image starting at the 
vertical center further illustrates differences. 
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The same spectral representation was then applied to the standard deviation filtered 
images, except the speckle contrast MTF served to reduce each component. Rearranging Eq. 
(1), 



° =K MTF {\fM I L = ZZ*^(|/~IK / L e ' 



(2nf?x+2wf?y) 



(6) 



The predicted speckle contrast image is then just the ratio: K pred = (Ti mag J(I image ). 
4. Results 

Figure 2(a) shows a human radial artery illuminated under different sections of a Gaussian 
beam. Raw images and generated flow maps are shown in Figs. 2(b) and 2(c). 
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Fig. 3. (a) Intensity image for three distinct beam shapes, a Gaussian shape, output through a 
frosted glass diffuser, and flat top. (b) Experimental speckle contrast images computed using 
7x7 neighborhoods, (c) Predicted speckle contrast using the forward cSFDI model based on the 
spatial frequencies of the beam shape. Note that artifacts within the beam shape are amplified 
due to the high pass nature of the speckle MTF. (d) (e) Effects from systematic variables such 
as vignetting, sliding window filter response, or low counts are shown to be negligible by 
direct comparison with a static object, where the speckle contrast is seen to be constant and 
approximately equal to /?, as expected. 
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One observes a bulk increase in calculated flow index of approximately 15% as the vessel 
branch is moved radially outward, away from the peak of the beam profile, as the photons 
exiting the tissue here have travelled longer average path lengths. This occurs because the 
measurement point is made further from the illumination peak at the center of the beam. 

This effect was studied more closely using homogenous liquid phantoms, where DFT's 
could be used to generate a forward contrast image using Monte Carlo simulations of cSFDI. 
Results are showed in Fig. 3 for three common beam profiles, a standard Gaussian shape, 
output through a diffusing glass plate, and a flat homogenous profile. The reflectance 
averaged over 500 exposures is shown in Fig. 3(a). Figure 3(b) show the experimental speckle 
contrast calculated with the usual sliding window standard deviation and mean filters. These 
are to be compared with the predicted values using cSFDI (Eq. (5) and Eq. (6)) shown in Fig. 
3(c). 

It is clear from the images that cSFDI is able to reasonably predict the behavior of the 
speckle contrast. The root mean square error (RMSE) was calculated at each pixel and 




Fig. 4. Images of mouse brain with a Gaussian beam shape (left) and flattened beam shape 
(right), (a) Raw images are shown where the shape of the beam is evident, (b) Flow index maps 
and (c) a specific ROI of super saggital sinus vessel illustrate a difference in perceived venous 
flow that depends on the beam shape, (d) The percentage difference between flow maps show 
significant variation across the field. 
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normalized by the theoretical value predicted by cSFDI to obtain a percent error. The average 
percent error over all pixels was less than 5% for all three shapes. 

To act as a control, speckle images were also taken of static scattering objects using a 
Spectralon reflectance standard. The raw images and calculated contrast are shown in Figs. 
3(d) and 3(e). As expected, the calculated contrast is effectively uniform across the image 
space regardless of beam shape and approximately equal to /?. 

Figure 4(a) shows the raw speckle image of a mouse brain for two intensity profiles, a 
Gaussian and flat shape. Speckle flow maps were generated and are shown in Fig. 4(b). 
Although the large variations within the image between vessels and tissue conceal the effect, 
inspection of a single section of the superior sagittal sinus vessel in Fig. 4(c) shows a bulk 
change of approximately 10% in calculated flow index between the two images. Finally, the 
percent difference image between the two flow maps is calculated across the cortical surface 
and shown in Fig. 4(d). 

Here, the differences in flow caused by the beam shape are clear, showing changes in 
some parts by as much as -25% near the peak of the Gaussian, due to shorter average path 
lengths, and + 25% looking further out along the tail. This agrees with predictions made by 
cSFDI and the in vivo human artery experiment presented above. 

5. Discussion and conclusion 

Our results show that the source intensity profile can have a significant effect on the speckle 
contrast and thus on the resulting speckle flow index. This phenomenon occurs because of 
variations in the average path length of photons across the beam profile. When measuring 
near a local intensity maximum, photons generally need to travel shorter path lengths to reach 
the detector, and thus transfer less overall momentum, as compared to lower intensity regions. 
The effect causes changes in the measured flow index across the field of view, and if not 
accounted for or minimized, would confound any attempt to determine absolute rates of blood 
flow. 

Under controlled conditions, the close agreement between the forward model predicted by 
cSFDI and experimental results suggest this technique is an accurate method for predicting 
the effect of light structure on speckle contrast. The relative ease with which the speckle 
contrast filter K MTF ij) can be applied to any complicated structure using DFT's, such as the 
"diffuser shape" in Fig. 3, illustrate the utility of the method in an imaging geometry. In 
principle cSFDI could be used to generate a multiplicative "correction" factor to account for 
beam shape. However, this creates a number of extra complications when imaging a 
heterogeneous sample. First, diffuse beam shape can no longer be analyzed using DFTs, as 
heterogeneities would affect the Fourier spectrum. In addition, the optical properties must be 
known at each location in the image and the speckle MTF and forward model must generated 
for each set. If the beam shape and bulk optical properties were assumed a priori, a simple 
correction factor could be generated, but this approach would contain errors that scale with 
the magnitude of the heterogeneities. Perhaps the ideal solution is to remove structure from 
the beam shape before measurements in the interest of quality data collection. Applying this 
philosophy to an in vivo experiment in a mouse brain shows major differences when taking 
steps to appropriately homogenize the beam. 

There are several methods to help reduce variation in the beam profile. The simplest is to 
expand the beam to a larger radius, effectively increasing the width of the natural Gaussian 
profile. The disadvantage to this method is a significant decrease in the intensity per unit area, 
and that the addition of optical elements can create more artifacts and aberrations. Another 
possibility is expanding the beam radius before entering a diffuser plate and decreasing the 
distance to the object. It is known that coherent light through a diffuser will project an 
intensity speckle pattern on an object some distance away, where the speckle size is governed 
by the incident spot radius and the projection distance. If this secondary intensity speckle 
pattern is fine enough such that the Nyquist criterion is not met, the path length distribution 
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will be uniform across the pixel array due to under sampling. For a circular spot of diameter D 
incident on a diffuser and projected a distance z away, the power spectral density is defined 
[26] to be a bandwidth limited function with maximum frequency f b : 

/»=-£-■ (7) 

AZ 

Thus, to avoid sampling the intensity speckle pattern with the camera array, the pixel 
width x p must be larger than the Nyquist limit: 

x> = . (8) 

' 2f b 2D 

Note this under sampling serves to blur the projected secondary speckle in conjunction with 
standard light diffusion through tissue, creating a very smooth profile. However, while the 
Gaussian profile of the beam is broadened when passing through a diffuser, it is not 
eliminated. 

A third option is to use a flat top beam shaper with a homogenizing element. Commercial 
ones are typically manufactured using a fly's eye micro-lens array, and these optical systems 
are also used in commercial projectors to form a homogenous rectangular beam. Differences 
in path length between elements of the lens array can cause a secondary interference speckle- 
like pattern [32], but our experience has been that these have been under sampled and blurred, 
as in Eq. (8), to create the homogenous flat shape seen in Fig. 3. 

In conclusion, the shape of the source beam causes a significant effect on the speckle 
contrast and speckle flow index measured in LSI. We demonstrated this effect using two types 
of in vivo measurements: small animal brain, and human radial artery. Using controlled 
microsphere phantoms, the effect was quantified and predicted using cSFDI framework. This 
was done both for a Gaussian beam and for the more complicated intensity profile from a 
diffusing glass plate. We showed that spatial frequency domain analysis provides a means to 
predict the effect of any beam shape on the speckle contrast by calculating a speckle contrast 
MTF. We also investigated and reported methods for increasing the beam homogeneity in LSI 
experiments. 
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